function cMap = makeColorMap(varargin)
%% MAKECOLORMAP makes smoothly varying colormaps
% a = makeColorMap(beginColor, middleColor, endColor, numSteps);
% a = makeColorMap(beginColor, endColor, numSteps);
% a = makeColorMap(beginColor, middleColor, endColor);
% a = makeColorMap(beginColor, endColor);
%
% all colors are specified as RGB triples
% numSteps is a scalar saying howmany points are in the colormap
%
% Examples:
%
% peaks;
% a = makeColorMap([1 0 0],[1 1 1],[0 0 1],40);
% colormap(a)
% colorbar
% 
% peaks;
% a = makeColorMap([1 0 0],[0 0 1],40);
% colormap(a)
% colorbar
% 
% peaks;
% a = makeColorMap([1 0 0],[1 1 1],[0 0 1]);
% colormap(a)
% colorbar
% 
% peaks;
% a = makeColorMap([1 0 0],[0 0 1]);
% colormap(a)
% colorbar

% Reference:
% A. Light & P.J. Bartlein, "The End of the Rainbow? Color Schemes for
% Improved Data Graphics," Eos,Vol. 85, No. 40, 5 October 2004.
% http://geography.uoregon.edu/datagraphics/EOS/Light&Bartlein_EOS2004.pdf

defaultNum = 128;
errorMessage = 'See help MAKECOLORMAP for correct input arguments';


sz = cellfun('prodofsize',varargin);

constraits = vertcat(varargin{sz == 3});
steps = [varargin{sz == 1}];

numcstr = size(constraits,1);
numsteps = numel(steps);

if numsteps < 1
  steps = round(diff(linspace(1,defaultNum+1,numcstr)));
elseif numsteps == 1
  steps = round(diff(linspace(1,steps+1,numcstr)));
elseif numsteps > numcstr-1
  steps = steps(1:numcstr-1);
elseif numcstr-1 ~= numsteps
  steps = round(diff(linspace(1,defaultNum+1,numcstr)));
end

steps(1:end-1) = steps(1:end-1)+1;

cMap = [];
for k=1:numcstr-1
 cMap = [cMap(1:end-1,:) ;
   interpMap(constraits(k,:), constraits(k+1,:), steps(k))];
end

% size(cMap)
% return
% if nargin == 2 %endPoints of colormap only
%     color.start  = varargin{1};
%     color.middle = [];
%     color.end    = varargin{2};
%     color.num    = defaultNum;
% elseif nargin == 4 %endPoints, midPoint, and N defined
%     color.start  = varargin{1};
%     color.middle = varargin{2};
%     color.end    = varargin{3};
%     color.num    = varargin{4};
% elseif nargin == 3 %endPoints and num OR endpoints and Mid
%     if numel(varargin{3}) == 3 %color
%         color.start  = varargin{1};
%         color.middle = varargin{2};
%         color.end    = varargin{3};
%         color.num    = defaultNum;
%     elseif numel(varargin{3}) == 1 %numPoints
%         color.start  = varargin{1};
%         color.middle = [];
%         color.end    = varargin{2};
%         color.num    = varargin{3};
%     else
%         error(errorMessage)
%     end
% else
%     error(errorMessage)
% end
%    
% if color.num <= 1
%     error(errorMessage)
% end
% 
% if isempty(color.middle) %no midPoint
%     cMap = interpMap(color.start, color.end, color.num);
% else %midpointDefined
%     [topN, botN] = sizePartialMaps(color.num);
%     cMapTop = interpMap(color.start, color.middle, topN);
%     cMapBot = interpMap(color.middle, color.end, botN);
%     cMap = [cMapTop(1:end-1,:); cMapBot];
% end
    

function cMap = interpMap(colorStart, colorEnd, n)

for i = 1:3
    cMap(1:n,i) = linspace(colorStart(i), colorEnd(i), n);
end

function [topN, botN] = sizePartialMaps(n)
n = n + 1;

topN =  ceil(n/2);
botN = floor(n/2);
% Copyright 2008 - 2009 The MathWorks, Inc.
